#pragma once

#include <tuple>
#include <math.h>

#include "Transformation.hpp"


// > 总岁差(general precession) => 由赤道岁差和黄道岁差组成
// ---------------------------------------------------------------------------------------------------------------------------------------------------------
// > 岁差(axial precession)
//   -> 黄道岁差(PrecessionInEcliptic)
//         - 描述：太阳系内其它行星质点引力对地球轨道面产生的摄动，使得黄道面的空间取向发生变化，引起春分点沿赤道运动;
//         - 符号：λ(lamdba),
//         - 起点：春分点
//         - 单位：度，更小的单位是“角分”和“角秒”
//         - 范围：[0°，360°] 由西向东量度
//   -> 赤道岁差(PrecessionInEquatorial) 
//         - 描述：由于日月和行星对地球赤道突起部分的较差力矩而引起的赤道面（极轴） 的长期进动;
//            -  可用近代的天文和测地技术(如VLBI、LLR等)的观测确定;
// 
//         - 符号：β(bata),
//         - 起点：黄经圈
//         - 单位：度，更小的单位是“角分”和“角秒”
//         - 范围：[-90°, +90°]北黄极为+90°，南南极为-90°(向北为正，向南为负)
// ---------------------------------------------------------------------------------------------------------------------------------------------------------

//T:为参考历元丁2000.0(丁D=2451545.0)至观测历元t之间的儒略世纪数：T = (JD(t) - 2451545.0) /36525.0

struct PrecessionTable
{
    // 赤道岁差(日，月岁差);
    double Psi[6];   // Ψ;(春分点在黄道上的向西移动，观测历元t的春分点移动量);//fi;
    double Omega[6]; // ω;                                                    //w;

    // 黄道岁差(行星岁差);
    double Pa[6];        // P;
    double Qa[6];        // Q;
    double Epsilon[6];   // ε: 为平黄赤交角变化量计算参数，计算公式如下(单位角秒);
    double Lambda[6];    // λ; 春分点在天球赤道上的东移量计算参数;

    //黄道岁差参数;
    double Eta[6];     // η;(自导出),根据P=sin(η)*sin(Π):起始历元黄道面到终历元黄道面的夹角;
    double Pi[6];      // Π;(自导出),根据Q=sin(η)*cos(Π)
    double p[6];       // ι;黄经总岁差(p表示)变化量计算参数;//E = E0 + E1*T + E2*T2 + E3*T3;

    // 赤道岁差参数;
    double Phi[6];       // φ;//X轴;            //th;theta;
    double Sigma[6];     // Σ:地球转动角;//Y轴; //Z;Zeta;
    double Zeta[6];      // ζ;//Z轴;            //z;zeta;
};

// IAU 岁差模型;
enum class PrecessionType: char
{
    IAU1978,
    IAU2000,
    IAUP03,
};
//岁差;
class Precession
{
private:
    static const PrecessionTable PrecessionTable1976;
    static const PrecessionTable PrecessionTable2000;
    static const PrecessionTable PrecessionTableP03;

    ~Precession() = delete;
    Precession() = delete;

    /**
     * @brief 计算岁差项;
     * 
     * @param T          儒略世纪数;
     * @param entry      计算条目数据;
     * @param entrySize  计算条目数据的大小;
     * @return double 计算后的数据;   
     */
    static double CalcEntry(double T, const double entry[],const size_t entrySize) noexcept;

    /**
     * @brief 获取岁差表;
     * 
     * @param type @PrecessionType 默认使用 PrecessionType::IAUP03
     * @return std::tuple<const PrecessionTable*, size_t> <表数据,表数据大小>;
     */
    static std::tuple<const PrecessionTable*, size_t> GetPrecessionTablePtr(const PrecessionType type = PrecessionType::IAUP03) noexcept;

    /**
     * @brief 岁差赤道旋转函数-<公用函数>;
     * 
     * @param llr   球面坐标;
     * @param phi   φ 值;
     * @param sigma Σ 值;
     * @param zeta  ζ 值;
     * @return Coordinate3D  旋转后的坐标;
     */
    static Coordinate3D CalcEquatorConversion(const Coordinate3D &llr, double phi, double sigma, double zeta) noexcept;
    // 没使用;
    /**
     * @brief 岁差黄道旋转函数-<公用函数>
     * 
     * @param llr   球面坐标;
     * @param eta   η 值;
     * @param pi    Π 值;
     * @param p     p 值;
     * @return Coordinate3D  旋转后的坐标;
     */
    static Coordinate3D CalcEclipticConversion(const Coordinate3D &llr, double eta, double pi, double p) noexcept;
public:

    //岁差旋转
    //J2000赤道转Date赤道
    /**
    * @brief 岁差旋转 J2000赤道球面坐标点 转 日期赤道球面坐标点
    * //CDllr_J2D
    * @param T
    *            儒略世纪数;
    * @param llr
    *            当前球面坐标;
    * @param type
    *            岁差类型 @see PrecessionType；
    *
    * @return 计算后的球面坐标;

    */
    static Coordinate3D ToDateEquator(const double T, const Coordinate3D llr, const PrecessionType type = PrecessionType::IAUP03) noexcept;

    //J2000赤道转Date赤道
    /**
    * @brief 岁差旋转 日期赤道球面坐标点 转 J2000赤道球面坐标点
    * //CDllr_D2J
    * @param T
    *            儒略世纪数;
    * @param llr
    *            当前球面坐标;
    * @param type
    *            岁差类型 @see PrecessionType；
    *
    * @return 计算后的球面坐标;

    */
    static Coordinate3D ToJ2000Equator(const double T, const Coordinate3D llr, const PrecessionType type = PrecessionType::IAUP03) noexcept;

    /**
    * @brief 岁差旋转 J2000黄道球面坐标点 转 日期黄道球面坐标点
    * //HDllr_J2D
    * @param T
    *            儒略世纪数;
    * @param llr
    *            当前球面坐标;
    * @param type
    *            岁差类型 @see PrecessionType；
    * 
    * @return 计算后的球面坐标;

    */
    static Coordinate3D  ToDateEcliptic(const double T, const Coordinate3D llr, const PrecessionType type = PrecessionType::IAUP03) noexcept;
  
    /**
    * @brief 岁差旋转 日期黄道球面坐标点 转 J2000黄道球面坐标点
    * //HDllr_D2J
    * @param T
    *            儒略世纪数;
    * @param llr
    *            当前球面坐标;
    * @param type
    *            岁差类型 @see PrecessionType；
    *
    * @return 计算后的球面坐标;
    */
    static Coordinate3D  ToJ2000Ecliptic(const double T, const Coordinate3D llr, const PrecessionType type = PrecessionType::IAUP03) noexcept;


    /**
     * @brief 计算黄赤交角;
     * [hcjj]
     * @param t 
     *          世纪数;
     * @return double P03黄赤交角;
     */
    static double EclipticObliquity(double t)
    {
        const double t2=t*t, t3=t2*t,t4=t3*t,t5=t4*t;
        return Transformation::ArcsecToRadians( (84381.4060 -46.836769*t -0.0001831*t2 +0.00200340*t3 -5.76e-7*t4 -4.34e-8*t5));
    }
};
